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Abstract 
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to estimate multigluon production processes at high energy hadron colliders. The 
partonic reactions gg — > (n — 2) g with n < 9 are studied and comparisons with 
known approximations are presented. 
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1 Introduction 



High-energy hadron colliders will be valuable tools for developing particle physics in the 
next decade, especially with the planned run of the Large Hadron Collider (LHC), which 
will operate at \/s = 14 TeV, the highest energy achieved so far. Multijet final states are 
of great importance at high energies, and in fact events with very high jet multiplicities, 
up to six jets, have already been observed and measured at the TEVATRON collider 

On the other hand, the estimate of the multijet production cross sections as well 
as their characteristic distributions is a formidable task from the theoretical point of 
view ||. Perturbation theory based on Feynman graphs, the 'brute-force' method, runs 
into computational problems, since the number of graphs contributing to the n-gluon 
amplitude grows like n!, jeopardising the computational efficiency In order to simplify 
the computation, colour decomposition [§, f|, [| has been used extensively in the past. 
According to this scheme, the n-gluon amplitude can be written as: 

M=2tg n - 2 ]T Tr{t ai ...t an )C{l,...,n) (1) 

P(2,...,n) 

where t a are the generators of the SU(N) group in the fundamental representation, g 
is the gauge coupling constant and C(l, . . . , n) are gauge invariant subamplitudes that 
do not contain any explicit information on the colour structure. The computational 
complexity is now due to the large number of gauge invariant subamplitudes, which is 
roughly proprtional to n\. Moreover the colour summation leads to a colour matrix with 
dimensionality (n — 2)! x (n — 2)! |J, another source of computational inefficiencies. 

On the other hand recursive methods [0, H, |[ || [K| have been used extensively in order 
to overcome the computational obstacles of the 'brute-force' method. Following this line 
of thinking we present in this paper results on gg — > (n — 2) g with n up to n = 9. The 
main ingredients we have used, in order to make the computation as simple and efficient 
as possible, are 

1. The amplitude is computed using recursive equations which result in a computa- 
tional cost growing asymptotically as 3™, to be compared with the n! growth in the 
'brute-force' method ]TT] . 

2. Colour and helicity structure are appropriately transformed so that a replacement 
of the usual summation with Monte Carlo integration is made possible. 



2 The Recursive Equations for QCD 
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2.1 The Generating Function 

Let us assume that (f>(x) represents a field with arbitrary quantum numbers. The gener- 
ating function of all connected Green functions is given by the standard formula [12 

Z[J] = -T[(j>] + J d 4 x(f)(x)J(x) , (2) 

subject to the condition that is a solution of the field equations 

T[0] is the effective action which at tree order is given by 

r[0] = J d 4 xC, (4) 

where £ is the Lagrangian of the theory. For a process with n external particles the 
source, J(x), is given as an expansion over the free asymptotic states, 

n 

J(x) = 5>e-^, (5) 
i=i 

or in momentum space by 

n 

J(?) = (27r)*£a,tf(g-p,), (6) 
i=i 

The standard LSZ reduction formula relates the amplitude with the n-th derivative of the 
generating function 

n cm y 

A(j h ,p 2 , ... )Pn ) = limnt^fa)]- 1 ^ dan (7) 
where [Z^ip)]' 1 is the inverse propagator and the on-shell limit is understood. Now it 



is not difficult to prove that the field equations aquire the following solution |T0| , |13 

2 n — 2 

0(g) = (2vr) 4 Y, b m 5(q-P m ) (8) 

m=l 

where P are all possible sums of the n external momenta, Pi, i — 1, . . . , n. For example, 
Pi + P2-, Pi + P4 + Pni etc. Their number is simply determined by the binomial sum 

n-l / \ 

n 



E( )=2 n - 2 



2 



The coefficients b m satisfy non-linear equations which can be solved iteratively leading 
to recursive equations. Since we are interested in multigluon amplitudes, we study the 
derivation of these equations in the case of QCD. The Lagrangian for QCD is given by: 

C = C + C int (9) 
Co = -\{d,Al - drAWPA" - PA*") - ^ A l? ( 10 ) 



,2 



Crnt = - 9 -f abC {d^A a v - d^A^A™ - Lf^efcde^^^ (n) 

where £ is the gauge fixing parameter. We shall take the value £ = 1 for the rest of this 
paper. In momentum space the effective action, r = / d 4 xC, is written as: 

1 1 q A a (-q)5 ab g^q 2 Al( q ) (12) 



(2tt) 4 " 

.2 ,. / 4 A<±„ \ I 



ff / (ft T^ll^Kft,..., q*)A afl (<li)A b »(q 2 )A™(q 3 ) ^(g 4 )(27r) 4 5 4 Q> 



8=1 



where V^(qi, q 2 , q-i) and G^^ x (qi, q^-, <?3, qi) are the usual three- and four-point vertices 
in QCD 0. 

For n-gluon scattering the source is given by 

n 

J a , = (2vr) 4 £ ^(Pi, \W ai 5{q - ju) (13) 

where {pf , Aj, a{\ are the momenta, helicities and colours of the external gluons. All 
momenta are taken to be incoming. The field equations aquire a solution of the form 

2"— 2 

^(g) = (27r) 4 £ b«(m)5(q - P m ) (14) 

m=l 

where P m are as in Eq.(|8|). The form of the solution, Eq.flT4]), suggests that a natural 
ordering of the momenta can be done in the following way: for the number m, lying 
between 1 and 2" - 2 we take its binary representation, which is an n-dimensional vector, 
{mi}, % = 1, . . . , n, with entries either or 1, and we define 

n 

PL = E m ^ ■ 

1=1 
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In the usual perturbative sense, one performs an expansion in terms of the coupling 
constant 

n— 1 

fc=i 

and the field equations take the form: 

& £K k) = ^ E Sm\ mim2 s k \k lk2 v;;^ (p m , p mi , p m2 )b% ( mi , kjv* (m 2 , h) 

E 6 m \rn im2mA k\k lk2 k A G a ^^\ (mi, fci)^ 2 2 (m 2 , fc 2 )i£ 3 3 (m 3 , fc 3 ) , 



6P 2 



(16) 



where summation over repeated indices is implied and 

The initial values for the iteration are given by the source term and read 

h a, -n _ / ^(P m ,X m )S aam m = 2 i -\i = l,...,n 
^ ' } ~\Q otherwise [U) 

The amplitude is then given by 

A = Y<b a {lA)-b a {2 n -2,n-l) (18) 

a 

where the "means that the propagator factor 1/P|„_ 2 has been removed. 
One can further simplify Eq.(|T6|) by defining a Lorentz scalar as 

b a (m,k,\) = e^P m ,\)b;(m,k) (19) 

where A stands for the helicity of the particle and takes the values A = ±, L, for transver- 
sally and longitudinally polarized gluons. Eq. ([[£]) can be inverted with the solution 

b«(m, fc) = E e»(P m , X)b a (m, k, A) + p J**^ £ ' Pm . (20) 

Ward identities guarantee that all terms proportional to P m drop out in the calculation 
of the physical amplitude, so that Eq. (|I6| ) can be written as 



b a (m,k,\) = — — ^ S m \ mim2 5k\k 1 k 2 V aaia2 (P rn , A; , P mi , Ai; P m2 , A 2 ) 
6 ai K^i,Ai)6 a2 (m 2 ,A; 2 ,A 2 ) 

c r> 2 ^ / ( ^m|mim2m3 ( ^fc|fciA;2fc3 
^ m 

6 ai (mi, fci, Ai)6 a2 (m 2 , k 2 , A 2 )6 a3 (m 3 , fc 3 , A 3 ) (21) 
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where 



and 



V abc (P t , A 4 ; P,, Xj] P k , A fc ) = e»(Pi, X l )e^P v X^(P k , \ k )V$&P u P,, P k ) (22) 



G abcd (X u A 2 , A3, A 4 ) = e"(P ls AxK(P 2 , A 2 )e K (P 3 , A 3 )e'(P 4 , A 4 )G^ P (23) 

Experience with scalar theory has shown that the computational cost grows like 3 n 
for a 3 theory and 4 n for a 4 theory [n]. This suggests that the computational cost 



will be lowered if we succeed in writing all equations in terms of three- vertex interactions. 
This can be done by introducing an auxiliary field X and replacing the four-gluon part 
of the Lagrangian with 

C x = -gf abc AlA b v X^ u c - X^X^ a . (24) 
Then Eq . (|2l|) is written as 

b a (m,k,X) = S m \ mim2 5 k \ klk2 V aaia2 (P m , A; , P mi , Ai; P m2 , A 2 ) 

6° 1 (mi,A;i,Ai)6 a3 (m2,A;2,A2) 

— ^ \ "* v r H^ U (P X- P X f aa i a 2 

/ y u m|mi77i2 "fc|fciA;2 JJ \ J mi / >; J mi j 

m 

b ai (m u k ll X 1 )X^{m 2l k 2 ) 
X^ u (m,k) = — - S m \ mim9 5 k \ k - lk9 H llu (P m - l , Ai; P m ^, A 2 )/ actl 



<J2 

/ I ' . . :,!'/. /•,/.. / /lll/l / / \ I . / , , . . . /v> I / 

6 ai (^i^i,Ai)6 a2 (m 2 ,A; 2 ,A 2 ) (25) 



where 

H^ipx, Xi,p 2 , A 2 ) = (e M (pi, Ai)e„(p2, A 2 ) - e„(pi, Ai)e /i (p 2 , A 2 )) . (26) 



2.2 Colour and helicity treatment 

Eq. (p5]) can be used to compute the n-gluon amplitude for an arbitrary momentum, colour 
and helicity configuration of the external particles. In order to have an estimate of the 
production probability, one has to sum over all colour and helicity configurations. Summa- 
tion over colours is a delicate subject. If one performs the summation in a straightforward 
way then one has to consider something like 8 n configurations for the n-gluon scattering. 
In this section we show how this summation can be replaced by integration, which is then 
suitable for Monte Carlo computation. 

As a first step a simplification of the colour structure of Eq . fl25|) is possible by defining 
the following object 
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b AB = J2t a AB b a , A, B = 1,2,3 (27) 

a=l 

where all other indices have been temporarily suppressed. The new objects are of course 
traceless three by three matrices in colour space. The interesting property of this colour 
representation is that it leads to a 'diagonalization' of the colour structure of the three- 
gluon vertex. More specifically, the colour part of the three-gluon vertex is now given 
by 

b b ^ 

f 1 ^AB^CD^EF = ~^($AD <>CF &EB ~ &AF &CB &Ed) (28) 

This colour structure is related to the colour flow occuring in the real physical process, 
where gluons can be represented by quark-anti-quark states in colour space and their self- 
interaction, as given by Eq.([28|), reflects the fact that colour remains unchanged on an 
uninterrupted colour line. 

Accordingly Eq. (|25|) is now transformed to 



b AB (m, k,X) 2 ^""^ ^m\mim2^k\kik2^ r {fmi A, ; ^mii Aj ; Pyrin A2) 

(b AC (m u fei, Ai)6 CB (m 2 , k 2 , A 2 ) - b AC (m 2 , k 2 , A 2 )&cb("Ii, h, Ai)) 
1 

^ y ^171 1 171 1 TU2 ^fc | fc 1 &2 ^ \U> ^ PfTl 1 A, Pm\ 1 ^l) 

m 

{b AC {mi, fci, \i)X£ B (m 2 , k 2 ) - X^ c (m 2 , k 2 )b C B{m 1 , h, A x )) 

^As( m ! ^) = 7 ^ ^mlmim^fclfcife^^mi ; Ajj -P m ^ , A 2 ) 

(&ac(wi, fci, Ai)6 CB (m 2 , k 2 , A 2 ) - b AC (m 2 , k 2 , \ 2 )bcB{m 1 , fc x , A x )) 

(29) 

where 

y{Pmi A; , -Pmn Ai; P m2 , A 2 ) 

stands for the momentum part of the three-gluon vertex, defined in Eq . (|22"D . 

The next step is to consider a colour representation so that summation over colours 
can be replaced by integration. We seek therefore a colour vector r] a (z) satisfying the 
following property: 

"[dz}r] a {z)r]\z) = 5 ab (30) 



with some properly defined measure [dz\. This can be realized by standard group- 
theoretical constructions. In the case of SU(3), we have to consider the 5-dimensional 
sphere [14}] parametrized by three complex numbers (zi,z 2 ,z 3 ) subject to the constraint : 

| Zl | 2 + | z 2 | 2 + | z 3 | 2 = 1 (31) 
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A useful parametrization of these three complex coordinates is as follows : 

Zl = e i<Pi cos 6 (32) 

z 2 = e i<t>2 sin^ cos £ (33) 

z 3 = e ^ 3 sin 9 sin f (34) 
O<0i<27r ; < 6> < f ; < £ < f 

and the group invariant measure is simply given by 

/ (lldzidsAsctzi^-l) (35) 

J \i=l ) i=l 

or, in terms of the polar variables, 

JJ J d<f>A J 2 d9 y 2 d£cosflsin 3 #cos£sin£ (36) 



\i=l 

It is now a matter of algebra to show that the vectors 

3 

V a (z) = v/24 £ z*(t a ) ijZj a = l,...,8 (37) 

satisfy Eq.(|3H|). 

As far as our recursive equations are concerned their structure remains unaffected and 
the only things to be changed are the initial conditions : 



8 f 1 \ 

b AB (m, 1, A) = foa ( m > A ) r / a (^) = ^AA m V^6 ZmA^B - ) 
a=l \ d / 



(38) 



where as usual m = 2 4_1 , i = 1, . . . , n, X m is the helicity and z m are the new continuous 
colour coordinates of the m-th gluon. 

In the same way, summation over helicity configurations of the external gluons can 
be replaced by an integration over a phase variable. This is achieved by introducing the 
polarization vector: 

ej(p) = y/2 (e%»(p, +) + e -*e"(p, -)) . (39) 
Then by integrating over we obtain the sum over helicities, 

- [* d<j>e%(p)(e;(p)y = £ ^(p,A)(^(p,A))*. (40) 

7T JO t i 



A=± 
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3 Results and Discussion 



Using the recursive equations described so far, we study the processes gg —>■ (n — 2)g, 
for up to n = 9. The Monte Carlo integration techniques used to represent the colour 
and helicity summation are expected to have an influence on the variance of the squared 
matrix element over the extended phase space, which now includes besides the momenta, 
the continuous colour variables, Zj, i = 1, 2, 3, and the helicity phase 0. In order to study 
the variance we have used the following computational schemes: 

1. By method I we mean that colour and helicity are treated as described in the 
previous section. 

2. In method II the helicity summation is replaced by a discrete Monte Carlo: in 
each 'event' we randomly choose a helicity configuration and then multiply by the 
appropriate combinatorial factor counting the total number of non- vanishing helicity 
configurations. 

3. Finally in method III a full summation over helicity configurations has been per- 
formed. 

All three methods have been used in the calculation of 3 and 4 gluon production. Moreover 
an 'analytic' approximation, the so-called SPHEL approximation [|T5|, TB|, has been used in 
the computation. 

In order to avoid collinear and soft singularities the following cuts have been used: 

p Ti >60GeV, |^| < 2, %>40° (41) 

For each given multiplicity, i.e. for a given value of n, we use a different centre-of-mass 
energy, y§, since higher final-state gluon multiplicities exhibit higher energy thresholds, 
due to the imposed cuts. The computation within a given gluon multiplicity utilizes the 
same set of momenta configurations, so a positive correlation is expected. On the other 
hand, the variance can be used to estimate directly the computational speed. As far as the 
CPU time is concerned, methods I and II use the same amount of resources, whereas in 
method III we need 10, for n = 5, and 15, for n = 6, times more CPU time. Finally the 
running strong coupling constant with Nj = 5, Aq C d = 200 MeV and a scale Q = PT, ma x 
has been used. 

In table [l] the total partonic cross section and the corresponding variance estimate are 
presented. We see that method I is generally more efficient than method II since the 
variance is always better. Moreover, it is also better than method III if we take into 
account the mulitplicative factor 10 (n = 5) and 15 (n = 6) needed for the computation. 
Colour and helicity integration do have a contribution to the variance as this can be 
inferred from the comparison with the variance computed in the case of SPHEL, where only 
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process 
1 



II III SPHEL 



2^ 3 2.43 ±0.18 2.15 ±0.19 2.35 ±0.1 2.43 ± 0.06 

400 GeV 0.034 ±0.01 0.036 ± 0.01 0.012 ± 0.002 0.004 ± 0.0003 

2-> 4 0.30 ±0.02 0.24 ±0.03 0.31 ± 0.02 0.40 ± 0.01 

600 GeV (4±1)10~ 4 (8±3)10" 4 (2±0.5)10~ 4 (2±0.3)10~ 4 

2->5 (4.7 ±0.3) 10~ 2 (7.6 ± 0.2) 10- 2 

900 GeV (9 ± 2) lO" 6 (6 ± 2) lO" 6 

2->6 (1.0 ±0.1) KT 2 (1.9 ± 0.1) 10- 2 

1200 GeV (1.0 ±0.4) 10~ 6 (1.0 ± 0.3) lO" 6 

2^7 (1.5 ±0.1) 10- 3 (2.3 ±0.1) 10~ 3 

1500 GeV (2.1 ±0.9) 10- 8 (1.3 ± 0.3) lO" 8 

Table 1: Total partonic cross section (in nb) and the corresponding variance with their 
Monte Carlo errors. 



phase-space integration is involved. Of course variance reduction schemes are possible but 
their study goes beyond the scope of the present publication. 

As a further step we have also considered differential distributions, such as the max- 
imum and minimum gluon transverse momentum, pT,max, PT,min, and the maximum and 
minimum invariant two-(gluons)jet mass. In Fig.[l| the four differential distributions are 
presented for the process gg — > 6g, whereas in FigfJ the same distributions are presented 
for the reaction gg — > 7g. In all qualitative agreement with SPHEL has been found. 

Note, however, that our computational scheme is an exact tree-order calculation of the 
amplitude, where no leading colour approximation nor special helicities selection have 
been used. 

In summary a computational method based on recursive equations has been presented 
which enables one to compute the amplitudes for the processes gg — ► (n— 2)g for at least up 
to n = 9. The method can be easily extended to include all partonic subprocesses [17]. It 
can provide a reliable tree-order event-generator for multi-jet production at LHC energies. 
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Figure 1: Differential distributions for the process gg — > 6g in pb/GeV; (a) maximum 
transverse momentum, (b) minimum transverse momentum, (c) maximum two-jet invari- 
ant mass and (d) minimum two-jet invariant mass. Black circles correspond to the exact 
result whereas triangles represent the predictions of SPHEL. 
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Figure 2: Same as Fig.ffl for gg — > 7g. 
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